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Theory of spin-orbit couphng in bilayer graphene is presented. The electronic band structure 
of the AB bilayer in the presence of spin-orbit coupling and a transverse electric field is calculated 
from first-principles using the linearized augmented plane wave method implemented in the WIEN2A; 
code. The first-principles results around the K points are fitted to a tight-binding model. The main 
conclusion is that the spin-orbit effects in bilayer graphene derive essentially from the single-layer 
spin-orbit coupling which comes almost solely from the d orbitals. The intrinsic spin-orbit splitting 
(anticrossing) around the K points is about 24 jieV for the low-energy valence and conduction bands, 
which are closest to the Fermi level, similarly as in the single layer graphene. An applied transverse 
electric field breaks space inversion symmetry and leads to an extrinsic (also called Bychkov-Rashba) 
spin-orbit splitting. This splitting is usually linearly proportional to the electric field. The pecu- 
liarity of graphene bilayer is that the low-energy bands remain split by 24 /xeV independently of the 
applied external field. The electric field, instead, opens a semiconducting band gap separating these 
low-energy bands. The remaining two high-energy bands are spin-split in proportion to the electric 
field; the proportionality coefficient is given by the second intrinsic spin-orbit coupling, whose value 
is 20 i-LeV. All the band-structure effects and their spin splittings can be explained by our tight- 
binding model, in which the spin-orbit Hamiltonian is derived from symmetry considerations. The 
magnitudes of intra- and interlayer couplings — their values are similar to the single-layer graphene 
ones — are determined by fitting to first-principles results. 
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I. INTRODUCTION 

Spin-orbit coupling is the most important interaction 
affecting electronic spin transport in nonmagnetic mate- 
rials. The use of graphene in spintronicai*^ would require 
detailed knowledge of graphene's spin-orbit coupling ef- 
fects, as well as discovering ways of increasing and con- 
trolling them. Bilayer graphene has attracted wide at- 
tention since it has a tunable gap caused by a transverse 
electric field (gate)<^"— Although the gap seems to satu- 
rate at around 0.3 eV, the possibility to turn on and off 
the electronic transport makes bilayer graphene a suit- 
able material for (mainly analog) electronic applications, 
including potentially spintronics ones. 

The electronic band structure of bilayer graphene de- 
rives from that of single-layer graphene, with taking into 
account for interlayer coupling. The basic electronic 
structure is well understood)^— What is not known yet 
is how, and by how much the realistic spin-orbit interac- 
tion modifies the electronic spectrum, especially when 
gated. The spin-orbit effects are important not only 
for the fundamental electronic band structure and its 
topology, but also for understanding such phenomena as 
spin relaxation (see the recent spin injection experiments 
in Refs. 013), spin Hall effect^ magnetoanisotropy, or 
weak (anti)localization. Conventional charge transport 
in bilayer graphene has been studied in detailfiS but 
spin transportii or spin-orbit induced charge electronic 
transporti^ in graphene bilayer is only starting to be ex- 
plored. 

In this paper we argue that the spin-orbit coupling 
in bilayer graphene comes mainly from the intra-layer 



spin-orbit coupling, contrasting earlier studie o^'^'^ that 
predicted large interlayer effects, enhancing the spin- 
orbit spectral splittings by a decade as compared with 
single layers. A recent model investigation^^ has re- 
ported on spectral features of a bilayer in the presence of 
both intra and interlayer hopping spin-orbit parameters. 
The single-layer spin-orbit physics comes from (nomi- 
nally unoccupied) carbon atom dxz ± ^dyz orbital a^^'^'' 
which hybridize with the Pz ones. The d orbitals spin- 
orbit coupling opens a gap at the K points — the corners 
of the hexagonal Brillouin zone — of the value of about 
24 fj,eVJ^^ The ct-tt hybridization, on which most stud- 
ies have focused, determines a further spin-orbit split- 
ting in the presence of an external transverse electric 
fieldr^iiiii^ that can also arise from the substratci^ The 
a-TT hybridization comprises the Stark effect (shift of the 
Pz orbitals in the presence of an electric field) , the on-site 
coupling of the shifted Pz and s orbitals, and, finally, the 
spin-orbit splitting of the in-plane Px ± ipy orbitals. 

As is now common we call intrinsic the spin-orbit 
splitting — which amounts to spectral anticrossings at and 
around K, while preserving the double spin degeneracy — 
in the absence of a transverse electric field, and extrinsic 
(or Bychkov-Rashbap^) the spin splitting (lifting of the 
spin degeneracy) in the presence of such a field. The dis- 
tinction is rather sharp in materials with a space inver- 
sion symmetry (such as bilayer graphene), in which the 
electronic bands arc always doubly degenerate; applica- 
tion of an electric field breaks the space inversion sym- 
metry and spin-orbit coupling removes this degeneracy,— 
usually in proportion to the field. An additional source 
of extrinsic spin-orbit coupling are adatoms;^ which may 
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cause spin relaxation in graphene by creating patches of 
enhanced spin-orbit couphngj^i^i^i 

This paper reports on comprehensive first-principles 
as well as tight-binding investigations of the electronic 
structure of bilayer graphene, in the presence of spin- 
orbit coupling, and in the presence (and absence) of 
an external transverse electric field. Furthermore, the 
most generic spin-orbit Hamiltonian consistent with the 
K points symmetry is derived for bilayer graphene in an 
external electric field. The first-principles method we 
use is the linearized augmented plane wave technique 
with the generalized gradient approximation^^ for the 
exchange-correlation potential, as embedded within the 
Wien2fc packagci^ The extended tight-binding model is 
constructed as an effective single-orbital hopping model 
reduced from a multi-orbital tight-binding scheme, us- 
ing only atomic spin-orbit coupling. The tight-binding 
parameters are obtained by fitting to the first-principles 
results around the K points. In the first-principles calcu- 
lations we take 1.42 A for the intra-layer atomic distance, 
and 3.35 A for the intcrlayer distance. The vacuum layer 
is taken to be 20 A, well enough to uncouple the bilayer in 
the transverse direction in the three-dimensional periodic 
structure calculation. 

Our main conclusion is that the single-layer spin- 
orbit coupling determines quantitatively the spin-orbit 
induced anticrossings and spin splittings at and around 
K(K'). The intcrlayer coupling of the two graphene 
sheets in the AB stacked bilayer produces parabolic 
bands around the K points. Two bands remain close to 
the Fermi level. These low-energy bands, one conduction 
and one valence, cross at K. Spin-orbit coupling leads to 
anticrossing of the two bands with the value of 24 /J.eV, 
as in single-layer graphene. This splitting is due to the 
presence of d orbitals in the 7r-bands. Removing d (and 
higher) orbitals from our calculation, the gap is reduced 
to about 1 /icV; the typical value coming from the ct-tt 
hybridization. In a transverse electric field an orbital gap 
opens, separating the conduction and valence bands. In 
addition, spin-orbit coupling leads to spin splitting, re- 
moving the spin degeneracy at a given momentum. This 
extrinsic splitting is peculiar in bilayers, due to intcrlayer 
orbital effects. At the K points, the spin splitting is inde- 
pendent of the electric field (at typical field magnitudes) , 
with the value of 24 fieY, given by the intrinsic splitting. 
Away from K the extrinsic spin-orbit coupling begins 
to dominate, giving the splittings of roughly 10 fieV per 
field of IV/nm, increasing linearly with increasing field. 
The fine structure of the spin splittings away from K is 
well described by the intra-layer spin-orbit couplings for 
the low-energy conduction band. Quantitative fits to the 
low-energy valence band and high-energy bands require 
introducing also intcrlayer spin-orbit coupling parame- 
ters. This we do by deriving the most general spin-orbit 
Hamiltonian at K, which has 10 real parameters. By em- 
bedding this Hamiltonian with the tight-binding scheme 
and fitting to our first-principles data we find that the 
intcrlayer parameters are in magnitudes similar (about 



10 ficV) to the intra-layer ones. Being off-diagonal, their 
actual contribution to the spectrum close to K is greatly 
suppressed. 

This work is organized as follows. In section II we 
present the tight binding model including the discus- 
sion of the relevant d orbitals for a general iV-layer AB- 
stacked graphene with intra-layer spin-orbit coupling. 
This model is discussed in detail for bilayer graphene 
in section HI. Next, in section IV we present the first- 
principles results and the fits from the tight-binding 
model around the K points. We discuss the intrinsic 
and extrinsic spin-orbit splittings of the bands and the 
intcrlayer spin-orbit couplings. Appendix A constructs 
the most generic effective spin-orbit Hamiltonian for the 
graphene bilayer in an external electric field at the K(K') 
point and, as well, for an arbitrary momentum k. 



II. MODEL HAMILTONIAN 

A. Tight-binding Hamiltonian 

The electronic structure of 7r-bands of graphite and 
of A^-layers-graphene is usually described by a tight- 
binding approximation, often parametrizied according to 
the Slonczewski-Weiss-McClure (SWMcC) model^i^l-^ 
and expressed in terms of the 7r-band on-site orbital 
Bloch wave functions: 

*A. (k) = ^ ^ &^(^+^^. ) pf [r - (R + t A J] , 

T (1) 
vl/B, (k) = ^ e"^(^+*-. ) pf [r - (R + tB J] , 

labeled by quasi- momentum k counted from the F point, 
sublattice pseudospin A or B, and the layer index i, which 
runs form 1 to (the number of layers). Here, and 
stand for the positions of the 2A^ atoms in the A^- layer 
elementary cell (for the AB-stacked bilayer graphene the 
situation is depicted at Fig. [T]) and the summation over 
R goes over all Bravais lattice vectors. 

The intra- and intcrlayer hoppings between the (effec- 
tive) Pz orbitals of the neighboring atoms are given by a 
set of parameters 7, schematically shown in Fig. [TJ Pa- 
rameters 7o and 71 describe the nearest neighbor intra- 
layer and intcrlayer hoppings, while 73 and 74 are indi- 
rect hoppings between the layers. In addition, the pa- 
rameter A is introduced to handle the asymmetries in 
the energy shifts of the corresponding bonding and anti- 
bonding states due to 71. The role of these hopping pa- 
rameters in the band structure and the correspondence 
between the conventional tight-binding model and the 
SWMcC parametrization is given in Ref. [13] ■ 

Using the ordered on-site orbital Bloch basis "^^Ai (k) , 
^'Bi(k), ^'A2(k), ^'Balk) the spinless 7r-band structure 
of the AB-stacked bilayer graphene with lower and up- 
per layers placed in potential V and —V, respectively, is 
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FIG. 1: (Color online) Left: sketch of the AB stacked bilayer 
graphene. The filled circles represent the carbon atoms on the 
sublattice (pseudospin) A and B, respectively. Index 1 stands 
for the lower and index 2 for the upper carbon mono-layers. 
Right: the hopping parameters 70, 71, 73, and 74 accord- 
ing to the Slonczewski-Weiss-McClure (SWMcC) convention. 
Atoms Bi and B2 are connected by 74 (not shown). 



described by the effective 4x4- Hamiltonian: 
HTB(k) = 



/ A + y 7o/(k) 74/* (k) 71 

70/* (k) +V 73/(k) 74/* (k) I 

74/(k) 73r(k) ~V 7o/(k) ' ' 

V 71 74/(k) 70/* (k) A-F 



Here /(k) is the nearest- neighbor structural function of 
the graphene hexagonal lattice with the lattice constant 
a = 2.46A: 



/(k) = e v^'^" [1 + 2e-'—^- cos(f fc,)] , 



(3) 



which is accommodated to our on-site tight-binding basis 
([1]) and chosen coordinate system, see Fig. [T] 

With the spin degree of freedom s = {t, 1} the on-site 
Bloch basis doubles: 



^'A.(k)^^'A.,.(k) = *A.(k)(»|5) 
^-B.lk) ^ *B.,s(k) = ^'B.(k) ® \s) 



(4) 



and the dimension of the TB Hamiltonian increases to 
4iV X AN: 



FTB(k) ^FTB(k)' 



ntxti 

I |;)(;| 



(5) 



In what follows when using the on-site Bloch states for 
the K point momentum K = (|7r/a, O) we employ short- 
handed notation: 



The TT-bands on-site wave functions of few-layer 
graphene built solely on the Pz orbitals are not affected 
by the atom's core spin-orbit L • s term, since the Pz 
orbitals carry zero orbital momentum. Therefore from 
the microscopical point of view, coupling of the Pz or- 
bitals to other atomic orbitals is needed to describe 
spin-orbit effects. The minimum realistic model employs 
d± = dxz ± irfyz orbitals, and also s and p± = p^ ± ipy 
orbitals if an external electric field is applied The re- 
sulting multi-orbital tight-binding model can be reduced 
by the Lowdin transformation^^ to obtain an effective 
Hamiltonian for the 7r-bands at the K point )^ where the 
normalized effective orbitals take the form: 



1 



[pz{A,) + ijd+{B,)] , 



pf{B,) = \ [p,(B,) +i7d-(A0] . 



(7) 



(6) 



The numerical value of the orbital mixing parameter 7 
was estimatedii to be 0.09. For our group-theory based 
analyzes of the spin-orbit effects at the K(K') point we 
do not need the explicit form of the effective orbitals, 
the only information we need is that they transform as 
TT-statcs. Exactly this requirement, as was already re- 
marked by Slonczcwski^ implies that pz orbital centered 
at atom A; (B;) should be paired with c?_(. (d-) orbital 
at atom B^ (A^). However, as we will see later, the ap- 
pearance of d± orbitals in p°^ is important for qualitative 
understanding of band spin splittings at K(K'). 



B. Spin-orbit Hamiltonian 

The spin-orbit Hamiltonian at the K(K') = (±|7r/a, O) 
point in the presence of an external transverse electric 
field is derived in detail from the group-theory argu- 
ments in Appendix A. It possesses single-layer-like in- 
trinsic and extrinsic (Bychkov-Rashba) spin-orbit cou- 
plings, whose strengths are given by the four intrinsic 
All, AJi, A12, X'i2 and two extrinsic Ao ± 2Abr param- 
eters, respectively. Their physical meanings and impor- 
tance will be explained later. Additionally, the symmetry 
group at K(K') allows four interlayer spin-orbit param- 
eters Ai, A3 and A4 ± SX4, whose indices refer to the 
spin-orbit interlayer geometry analogous to the SWMcC 
hopping convention. Within the basis of the on-site spin 
Bloch functions ^ the spin-orbit at K(K') Hamiltonian 
is given as follows: 
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FIG. 2: Calculated band structure of the 7r-bands of bilayer graphene along the FKM lines, where K = |rK| = 47r/(3a) with 
a — 2.46 A. First-principles results are shown by circles while tight-binding calculations by solid lines, (a) Low and high-energy 
bands, (b) Fine structure of the low-energy bands, (c) Detail view of the low-energy bands showing spin-orbit coupling induced 
anti-crossings at the K point and at A; = — 0.063 nm~^. 
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-i(Ao + 2Abr)s; 
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In the above expression r = 1(— 1) for the K(K') point, 
respectively, and the matrices Sz and = \{sx ± iTSy) 
stand for the z-component and raising and lowering spin 
operators. With zero electric field the bilayer graphene 
has inversion symmetry and the number of spin-orbit A- 
parameters reduces from ten to four, see the discussion in 
Appendix A, or Ref. [ll]; these are An = Aj;^, A12 = AJ2, 
Ao and A4 and all others are forced to be zero from sym- 
metry requirements. The spin-orbit couplings and hence 
Hso a-i'c momentum independent. However, at the end 
of App. A we discuss their possible k-dcpendcnt exten- 
sion, which, as we will see, plays a very minor role for the 
spectra near the K(K') point. Therefore the full model 
Hamiltonian for gated bilayer graphene in the vicinity of 
the K(K') point is: 



i7eff(k) = ffTB(k) ® I ' ' I 1 1 " 1 1 + i/so . (9) 



The resulting electronic spectra of bilayer graphene de- 
rived from Hamiltonian ^ in the presence of a trans- 
verse electric field, as well as our first-principles results, 
are presented in the following sections. 



III. BILAYER GRAPHENE 

The electronic bands of bilayer graphene around the 
K point are parabolic, in contrast to (in the absence of 
spin-orbit coupling) linear bands in single-layer graphene. 
When gated, a tunable band gap opens<^i^i^ It has re- 
cently been proposed that the effects of spin-orbit cou- 
pling in bilayer graphene are of the order of hundreds of 
/ieV, caused by effective spin-dependent interlayer hop- 
ping between p orbitals J^ii^ Our results presented below 
do not support this view. 



A. Summary of results 

Figure [5] shows the calculated electronic band struc- 
ture of bilayer graphene around the K point along the 
FKM high-symmetry lines; spin-orbit coupling is taken 
into account, but there is no applied electric field so the 
spin degeneracy is present. There is excellent agreement 
between the tight-binding model and the first-principles 
calculations, in all the energy scales shown. Each scale 
in the figure has its own physics. Figure shows 
the usual picture of four, spin-degenerate, parabolic tt- 
bands. The high-energy conduction and valence bands 
arc formed mainly from the Pz orbitals at atoms Ai and 
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FIG. 3: Bilayer graphene essentials in an external electric field 
obtained from first-principles calculations (circles). The fig- 
ure plots the electric field dependencies of (a) the electrostatic 
potential 2V; the slope is described by the effective interlayer 
distance of dcfi = 0.1 nm, which is defined by 2V — eEdcs; 
the dashed line here is the linear fit. (b) the hopping pa- 
rameter 7i obtained by fitting to the first-principles data at 
K; (c) the energy gap in biased bilayer graphene, compared 
to the voltage 2V (dashed-dotted line); the solid line here 
is the tight-binding calculation using the potential 2V from 
the first-principles data in (a); (d) the parameter A obtained 
by fitting to the first-principles data at K. Note that E is 
the actual external electric field and not the screened one, as 
presented for example in Ref . 0] . 



B2 at the K point (see Fig. [T]) . These bands are shifted 
in energy by about 340 meV by the direct interlayer hop- 
ping 71 away from the low-energy bands, formed predom- 
inantly by the orbitals at atoms A2 and Bi. The two low- 
energy bands, again one conduction and one valence, are 
closest to the Fermi level. In the tight-binding model the 
difference in the energy shift between the conduction and 
valence bands is taken into account by the parameter A. 
Although we do not explicitly specify so, the spectra at 
the K' points are identical and our discussion is valid also 
for them. 

Figure HJb) reveals a fine structure of the low-energy 
bands. The bands form two overlapping parabolas, cross- 
ing at the K point, directly at the Fermi level, as well as 
at the point of accidental crossing along the FK line, at 
about 0.5 meV above the Fermi level. This is the man- 
ifestation of the trigonal warping, that induces a break- 
ing of the Fermi surface in the vicinity of each Dirac 
point into four pockets (Lifshitz transition), see for ex- 
ample HI]. These crossings are governed by the indirect 
interlayer hopping parameters 73 and 74, which pull the 
two bands together. The spin-orbit coupling causes anti- 
crossings of 24 /icV, just as in single-layer graphene^ 
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FIG. 4: The positions of the lattice momenta along the FK 
line (measured from K), corresponding to the minimum of the 
conduction band (circles) and to the maximum of the valence 
band (squares), for different values of the transverse electric 
field. The first-principles and tight-binding results coincides. 



as seen in Fig. Wis). The anti-crossings collapse below 
1 /icV if d and higher orbitals are excluded from the cal- 
culations. 

Applying a transverse external electric field i? to a bi- 
layer places the two layers at a different electrostatic po- 
tential. In the tight-binding model this is described by 
introducing a potential 2V ^ which includes all possible 
screening effects and corresponds to the splitting of the 
low-energy bands at the K point.— Figure [3l^a) shows 2V , 
extracted from the first-principles calculations, as a func- 
tion of the electric field. The dependence is almost lin- 
ear with the slope of about 0.1 electron nanometers, cor- 
responding to the effective electrostatic bilayer distance 
dgff « 0.1 nm. The electric field induces also a slight 
variation of the parameters A and of the direct inter- 
layer hopping 71 , obtained by fitting to the tight-binding 
model. The corresponding dependencies are shown in 
Figs 131b) and[3Id). In Fig. |31[c) we plot the spectral gap 
as a function of the electric field. At low electric fields the 
band gap is manifestly indirect, with a marked difference 
between the minimum of the conduction band and the 
maximum of the valence band. As the electric field in- 
creases beyond 1 V/nm, the gap becomes almost direct. 
The lattice momenta of the conduction band minimum 
and the valence band maximum are plotted in Fig. |4] for 
reference. The corresponding physics of the gap opening 
is discussed below. 

At small electric fields, less then 6mV/nm, bilayer 
graphene is a semimetal. A finite Fermi surface of a tri- 
angular shape is formed from the low-energy bands; the 
electric field induces only small energy gaps at the cross- 
ing points of the two overlapping parabolas (manifesta- 
tion of trigonal warping). A further increase of the elec- 
tric field opens an indirect band gap between the maxi- 
mum of the valence band present at the KF line and the 
minimum of the conduction band present at the K point 
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FIG. 5: Calculated band structure of the vr-bands of bilayer graphene in a transverse electric field of 25mV/nm. The spectra 
along the FKM line are shown, with the K point and the Fermi energy at the origin; K = |FK| = 47r/(3a) with a = 2.46 A. 
The first-principles calculations are shown by circles, while the lines come from tight-binding modeling, (a) Low and high- 
energy bands, (b) The fine structure of the low-energy bands with the corresponding indirect band gap of 0.7 meV between 
k = —0.068 nm~^ and the K point, (c) Detail view at the low-energy conduction band split by extrinsic spin-orbit coupling; 
the maximum value of the splitting of 2Ai = 24yueV is at the K point and at fc = —0.063 nm~^. 




FIG. 6: Calculated band structure of the 7r-bands of bilayer graphene with the applied electric field of 1 V/nm. Circles show the 
results of the first-principles and lines of the tight-binding calculations, (a) Low and high-energy bands, displaying a bandgap. 
(b) View of the low-energy bands showing the (mostly) direct bandgap of 94.5 meV between the valence (at k = —0.2 nm~^) 
and the conduction (at k — —0.15 nm~^) bands, (c) Detail view of the low-energy conduction band showing the spin splitting, 
with the maximum value of 2Ai = 24 ^eV at the K point and in its close vicinity. 



(Lifshitz transition). The corresponding band structure 
is shown in Fig. [5j The global picture of the bands is seen 
in Fig. [5ja), while the opening of the indirect band gap 
of 0.7 meV is seen in the fine structure zoom in Fig.[5l[b). 

At electric fields greater than 0.45 V/nm the states 
at the K point repel significantly; the minima and the 
maxima of the valence and conduction bands arc present 
along the FK and KM lines, and not at the K points^ii^ 
The band gap becomes (mostly) direct; the minimum of 
the conduction low-energy band and the maximum of the 
valence low-energy band are in general at very close but 
still different momenta. The spectrum for the electric 



field of 1 V/nm is shown in Fig.[ni[a), with the fine struc- 
ture showing the opening of the direct gap of 94.5 meV in 
Fig-Eljb). Finally, an extreme case, that of a very high- 
field of 10 V/nm, is shown in Fig. [71 The direct band gap 
of 265 mcV is seen in Fig. [Tljb). Due to the shift of the 
conduction band minimum away from the K point, the 
bandgap at fields greater than 0.45 V/nm is no longer 
proportional to the potential 2V (which determines the 
splitting at the K point), but rather it saturates to a 
value of about 265 meV— i^i^ as is shown in Fig. ^c) . 

The applied electric field breaks space inversion sym- 
metry and lifts the spin degeneracy. The spin splittings 
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FIG. 7: Calculated band structure of the vr-bands of bilayer graphene with the applied electric field of lOV/nm. Circles 
show the results of the first-principles and lines of the tight-binding calculations, (a) Hyperbolic high-energy bands and the 
Maxican-hat structure of the low-energy bands, (b) The fine structure of the low-energy bands in the bandgap region with 
the (mostly) direct bandgap of 265 meV; the maximum of the valence band is at = —0.910 nm^^ and the minimum of the 
conduction band is at fc = —0.904 nm~^. (c) Detail view at the low-energy conduction band shows its splitting due spin-orbit 
coupling with the value of 2Ai — 24 /leV at the K point. 



for the low-energy conduction bands round the K points 
for E = 25mV/nm, E = 1 V/nm, and E = lOV/nm are 
shown in Figs ^c) , ^c) , and ^c) , respectively. At these 
large field ranges, the spin splitting at the K point is seen 
to be independent of the field, having a constant value 
of the intrinsic splitting of 24 /xeV, deriving from the d 
orbitals. 

As seen from Figs [21 El [H and [71 all the spectra around 
the K points, including the fine structures and the spin 
splittings, can be faithfully described by tight-binding 
modeling. In the following we analyze the spin-orbit- 
coupling induced anticrossings and spin splittings. 



low, reads 



B. Analysis of results 



Spin-orbit coupling at and near the K point 



The spin-orbit effects at K can be understood essen- 
tially in terms of the interplay between the electrostatic 
potential 2V, direct interlayer hopping 71, and intrinsic 
spin-orbit couplings controlled by An, A^, A^, Ajj. Ba- 
sically with only those parameters, the energy spectrum 
of bilayer graphene at the K point ordered from high to 
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(15) 
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as sketched in Fig. [H The above spectrum can be derived 
from the Hamiltonian ^ when treating the spin-orbit 
interaction in the first order perturbation theory. 

The values for 2V, 71, A, intrinsic spin-orbit cou- 
plings An, AJi, A12, AJ2 and direct interlayer spin-orbit 
parameter Ai are obtained by comparing the eigenvalues 
Eq. pom7p of the effective tight-binding bilayer Hamilto- 
nian, Eq.([9|), to the first-principles spectra at the K point. 
This analysis shows that the spin splittings of the low- 
energy valence and conductance bands at the K point are 
the same and do not depend on the applied electric field; 
the spin-orbit parameters for these bands are predicted to 
be 2An — 2AJi = 24 fieV. In contrast, the spin splittings 
of the high-energy valence and conductance bands at the 
K point depend on the applied electric field. However, 
the high-energy intrinsic spin-orbit couplings governing 
these splittings are field independent and their values are 
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FIG. 8: Energy spectrum diagram of the AB stacked bilayer 
graphene at the K point. The diagram at left shows the spin- 
orbit energetics of two interlayer-non-interacting (71 = 0) 
graphene sheets placed in the transverse potential difference 
2V . The diagram at right represents the spin-orbit energet- 
ics in the presence of a transverse external field including 
the direct interlayer interaction (in TB Hamiltonian mediated 
via the parameter 71). The states formed predominantly by 
atoms A2 and Bi are split by the intrinsic spin-orbit cou- 
pling 2Aii. These states form the low-energy valence and 
conductance bands. The states residing mainly on Ai and 
B2 form the high-energy bands and are shifted in the energy 
spectrum by the direct interlayer hopping 71. They are spin 
split by 2A12/V1 + (7i/V)^ The energies el (i = {1,2,3,4} 
and s = {t, 4-}) and the corresponding eigenstates are ordered 
from top to bottom; see Eqs. (|10II17|) . 



fixed by 2Ai2 — 2Aj2 = 20 /icV. Finally the direct spin- 
dependent interlayer parameter Ai = 0. The remaining 
hopping (73 and 74) and spin-orbit (Aq, 2Abr, A3, A4 and 
(5A4) parameters, as discussed in the next subsection, are 
chosen to reproduce the band-structure in the vicinity of 
the K point (see the largest ranges shown in Figs. [21 [SJ 
[6l and[7|). The parameters for selected values of electric 
field used in this paper are listed in Tab. HI 

Within the first order perturbation theory at the K 
point the eigenstates of the Hamiltonian ^ can be ex- 
pressed (apart from the overall normalization) in terms 
of the on-site spin Bloch wave- functions (|4]), for the form 
of the eigenstates see Fig. El where 



sj-ll + V^±V 
71 



(18) 



To understand qualitatively the bilayer spectrum and 
the spin splittings at the K point we need to (i) ap- 
proximate the spin-orbit interaction j^^r^i^V x p) • s 
with the bilayer graphene point symmetry by the atomic 
(and hence isotropic) spin-orbit interaction ^ L • s and 



(ii) take into account the d± states in the effective p'^^ 
orbitals, which enter the on-site Bloch wave functions, 
see Eqs. ([1]), ([4]) and The on-site Bloch wave func- 
tions carrying the pseudospin A are formed by d+ or- 
bitals contrary the on-site wave functions labeled by the 
pseudospin B which are composed of (i_ orbitals. The 
atomic spin-orbit energies of d+® |t) and d-® are 
equal and higher by 2^^ than the energies of d+ig) \\^) 
and d-® W). Hence we expect the general spectral 
tendency — the spin-orbit pulling the on-site Bloch states 
^'A,t and ^B,4, higher in energy by 47^^^ compared to 
"^A.i and vj/g ^f.. 

Let us firstly analyze the zero external electric field 
case without the spin-orbit interaction. When two single- 
layer graphene sheets are brought together, the Bloch or- 
bital states on the interlayer-direct-contact atoms Ai and 
B2 start to interact (via spin independent parameter 71) 
and they form the high-energy antibonding (conduction) 
and bonding (valence) states ^'ai + and — ^'Bsi 
respectively. The low-energy states, on the other hand, 
are formed by the indirect-contact Bloch orbitals ^I^Bi 
and 5'a2- If we now turn on the spin-orbit interaction 
and count the spin degrees of freedom we would see the 
following changes in the bilayer energetics: The two anti- 
bonding (bonding) states ^'Ai,t + *B2,t ^nd 'I'Ai.i+^Ba,^ 
(^'Ai,t - *B2,t and 'I'Ai.j, - ^Ba.i) should stay spin un- 
split since the opposite pseudospin components ^a,s and 
^'b.s entering the antibonding (bonding) wave functions 
are shifted opposite in energy by ^ L • s and hence there 
is no net spin splitting. The situation is different for the 
four low-energy states ^Bi,s and ^'a2,s- Their four-fold 
degeneracy is partially lifted when the spin-orbit inter- 
action is turned on; states ^'bi,^ and ^'A2,t reniain de- 
generate and become shifted in energy higher than the 
degenerate pair ^'Bi.t and ^'a2,^- This splitting is seen 
in the spectrum as the K point anticrossing (see Fig. [2]), 
and according to the above qualitative model we can fix 
the value of 47^^^ to 24 /leV. This reasoning is fully 
consistent with the tight-binding energy spectrum when 
plugging in Eqs. (jlOlllTp zero for V and Ai and setting 
All = Aj]^. The spectral situation is schematically de- 
picted on Fig. [Sj 

An external electric field breaks inversion symmetry 
and causes external spin splittings we observe in the 
first-principles spectra. Since the first layer is placed to 
the potential V and the second layer to the potential 
—V, we separate the spin-spht states 4'bi,4, ^nd ^'Bi,t 
away from the spin-split states ^'A2.t ^nd 5'a2.4,j what 
manifests in the low-energy spectrum as the band gap 
opening, see Figs. [SJ71 The energetics of the high-energy 
bands is somewhat different. The applied electric field 
affects the antibonding (bonding) Bloch states by raising 
the relative magnitude of the pseudospin component ^'ai 
(^'b2) over the component ^'b2 (^aJ, see the behavior 
of a± given by Eq. (|18p . The corresponding spin split- 
ting is then dictated by the dominant pseudospin orbital 
^Ai (^B2) when coupled to the spin, i.e. 5'Ai,t (^62,^) 
goes in energy higher than ^'Ai,^ (^62. t)- This ten- 
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74 


2V 


2Aii 


2Ai2 


Ao 


2Abr 


Ai 


A3 


A4 


5A4 


SO 


SWMcC 


A — 72 + 75 


70 


71 


73 


-74/2 


2V 






















0.0097 


2.6 


0.339 


0.28 


-0.140 





24 


20 


5 











-12 







E = 20 m V / nm 


n c\c\(\Ci 

0.009o 


2.6 


0.339 


0.28 


-0.145 


0.0013 


24 


20 


5 


0.25 





0.038 


-12 


-0.075 




E = 1 V/nm 


0.0096 


2.6 


0.339 


0.25 


-0.165 


0.1059 


24 


20 


5 


10 





1.5 


-12 


-3 




E = 6 V/nm 


0.0094 


2.6 


0.343 


0.29 


-0.143 


0.6238 


24 


20 


5 


60 





9 


-12 


-18 




£ = 10 V/nm 


0.0092 


2.6 


0.348 


0.26 


-0.100 


0.9572 


24 


20 


5 


100 





15 


-12 


-30 





TABLE I: Tight-binding (TB) parameters in the units of eVs and spin-orbit (SO) couplings in the units of /^eVs, obtained by 
fitting the band structure to the first-principles calculations. The signs of the TB parameters are chosen to be consistent with 
the Slonczewski-Weiss-McClure (SWMcC) parametrization,— which is also shown. The translation table of the parameters in 
the tight-binding and SWMcC models is obtained from band-structure fitting of graphite. The presented values of the TB 
parameters are similar to those found elsewhere,— "^^^^i^ and are consistent with the values of Ref. obtained from bilayer 
band-structure calculation using the WIEN2fc code. 




FIG. 9: Calculated spin splittings of the low-energy- 
conduction band for the electric fields of _E = 25mV/nm, 
E = IV/nm, and E — 6 V/nm. The solid lines are first- 
principles and dashed lines tight-binding results with the 
intra-layer spin-orbit effects. The arrows and circles indicate 
the positions of the conduction band minima. 



dcncy is again well confirmed by our tight-binding model 
spectrum, Eqs. (|10m7|) . and its eigenstates distribution. 
Fig. m and as well by the first-principles calculations. 



2. Interlayer spin-orbit couplings 

The intra-layer spin-orbit couplings An, A12 and Abr 
suffice to explain the bilayer spectrum in the presence 
of a transverse electric field directly at the K points. 
Around K points, including the regions of the valence 
band maxima and the conduction band minima, the de- 
scription is satisfactory for the conduction band only, see 
Fig. ini However, if we look at the low-energy valence 
band (or high-energy bands), we see that the fine features 
of the spin splittings differ from what the tight-binding 
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FIG. 10: Effects of the intra-layer spin-orbit couplings. Cal- 
culated spin splittings of the low-energy valence (left) and 
conduction (right) bands. First-principles results are shown 
as solid lines, while the tight-binding fits with intra-layer 
spin-orbit parameters only are dashed. The maxima of the 
valence bands and the minima of the conduction bands are 
indicated by thin dashed vertical lines. The conduction band 
spin physics appears satisfactorily described, but the valence 
band splittings are rather off for the tight-binding model. 



model with intra-layer spin-orbit coupling predicts. This 
is markedly seen in Fig. (TU] 

In the following we include into the picture interlayer 
spin-orbit couplings, motivated by our symmetry-derived 
Hamiltonian in App. A, and demonstrate a very good 
quantitative agreement with first-principles data. We 
stress that, (i) the interlayer spin-orbit couplings are of 
the same order as the intra-layer ones, that is typically 
10/LteV, and (ii) our fitting, while physically motivated 
and robust, can be in principle non-unique as there are 
in principle 10 parameters entering the spin-orbit Hamil- 
tonian, see Eq. ([S]). As such, the presented model, de- 
scribed by Hamiltonian ([5]), should be considered as a 
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FIG. 11: Best fitting tight-binding model. Calculated spin 
splittings of the low-energy valence (left) and conduction 
(right) bands. First-principles results are shown as solid 
lines, while the tight-binding fits, including all possible k- 
independent spin-orbit parameters, are dashed. The maxima 
of the valence bands and the minima of the conduction bands 
are indicated by thin dashed vertical lines. 

physically reasonable convenient minimum quantitative 
description of the spin-orbit physics in bilayer graphene 
with broken space inversion symmetry by a transverse 
electric field. If one aims to describe the low-energy con- 
duction band only, one can neglect these interlayer spin- 
orbit couplings entirely. 

Wc now explain the fitting procedure. The spin-orbit 
Hamiltonian in the external electric field, as presented 
by Eq. ([5]), or Tab. HI] in App. A, contains 10 real A pa- 
rameters. The diagonal ones, the intrinsic Ai, determine 
the spin splittings at the K point. They can be easily 
fixed, and we find that they do not depend on the electric 
field significantly (within one fieV). We then look at the 
spectral changes of the spin splittings around K point as 
different off-diagonal parameters vary, and compare the 
trends in the changes with the actual first-principles data. 
Moreover, we assume that 2Abr,, Ai, A3 and 6X4^ which 
are absent in zero electric field, scale linearly with the 
intensity of the applied electric field. One can see that 
Ai « 0, as changes in this parameter distort the picture 
away from the first-principles results. The direct inter- 
layer coupling is then largely spin-independent (governed 
solely by 71). _ _ 

The parameters Aq = Aq + 2Abr and Aq = Ao — 2Abr 
describe both the global and local breaking of space in- 
version symmetry. In the absence of an electric field, 
Aq = Aq = Aq, since 2Abr vanishes, see App. A. The pa- 
rameter Ao describes a local bulk-inversion-asymmetry 
physics: the electrons in one layer feel an effective elec- 
tric field due to the presence of the other layer. This field 
gives rise to a "local Dresselhaus"— spin-orbit coupling. 
Naturally, the field is opposite in the two layers so the net 
effect is zero, as there is no global bulk inversion asym- 
metry. We can estimate Aq along the following lines. The 



FIG. 12: Calculated spin splittings of the high-energy valence 
(left) and conduction (right) bands. First-principles results 
are shown as solid lines, while the tight-binding calculations, 
using the best fitting parameters, are dashed. 

direct interlayer coupling is given by the energy 71 « 0.3 
eV. Since the distance between the two layers is about 0.3 
nm, the effective electric field felt by each layer due to the 
presence of the other is about 1 V/nm. We know from 
single-layer graphene that such a field gives the (real) 
Bychkov-Rashba splitting of 2Abr = 10 ^eV. This gives 
an order of magnitude estimate Aq ~ 10 fJ.eV, which gives 
also a check on how reasonable the actual fit is. 

In our fitting, knowing both the complexity and crude- 
ness of the procedure, wc focused on obtaining both 
reasonable and robust results. To that end we made 
a restricted least-squares fit to the first-principles spin- 
splittings using the data for the electric fields of 1 and 6 
V/nm. The fits were simultaneous to both data sets, for 
the valence and conduction low-energy bands only. The 
fitting was restricted to the spectrum around K within 
2.5% in the two directions K-F and K-M, constraining 
the fits to closely preserve the extremal points (minima 
and maxima) of the spectral splittings, so that the over- 
all shape was correct. The crucial test of the robustness 
of the obtained parameters, and of the assumption of the 
linearity (in E) of the spin-orbit parameters, was (a) re- 
producing the spin splittings of the high-energy bands 
which were not used in the fitting procedure, and (b), re- 
producing the spin splittings of the low- and high-energy 
bands at intermediate electric fields, which were also not 
included in the fitting. The actual parameters of the fits 
are given as follows: 

Ao = 5 /ieV , 2Abr = 10 x £;[V/nm] /zeV , 

Ai = , A3 = 1.5 X E[Y/nm\ fieY , (19) 

A4 = -12 ^eV , 5X4 = -3 X £;[V/nm] ^jcV , 

where the numerical value of the electric field intensity 
E should be taken in the units of V/nm. 

The first-principles data and the tight-binding fits 
for the low-energy bands at electric fields of 1, 4, and 



11 




(3) 



FIG. 13: (Color online) Calculated, from tight-binding, the 
spin orientation of the low-energy conduction-band states as 
function of the momentum k for different directions for the 
electric field of E = 1 V/nm. The axis corresponds to 
the TKM line. The circles corresponds to (1) |k|/K= 0.01, 
with the corresponding energy e = 56meV and the angle 
between the spin pointing vector and the kz axis 9 — 5.5°; 
(2) |k|/K= 0.0225, e = 119 meV, and 6 = 45.5°; and (3) 
lkl/K= 0.05, £ = 355 meV, and 9 = 80.5°. 



6 V/nm, are shown in Fig. [TT] The first-principles re- 
sults for 4 V/nm is, as discussed above, were not used in 
the fitting of spin-orbit couplings (calculations for other 
fields give similar level of agreement). The fact that the 
first-principles results are well reproduced signifies the 
robustness of the procedure and validity of our assump- 
tions. The spin splittings in bilayer graphcne are rather 
complex, also considering that the interesting points arc 
not really K but the positions of the valence band maxima 
and conduction bands minima. There is a clear competi- 
tion between the intrinsic splitting Ai, dominating close 
to K, and the extrinsic (off-diagonal) splittings, dominat- 
ing at momenta away from the K point. The spin pattern 
of the low energy conduction band is shown in Fig. [T51 
at the K point the spin quantization axis is along z. This 
is due to the intrinsic spin-orbit coupling. Away from 
the K point the spin quantization axis is in the plane, 
reflecting the dominance of extrinsic SOC. 

The splittings of the high-energy bands, which nor- 
mally do not play a role in transport, provide additional 
assurance in the fitting. The results from our first- 
principles and tight-binding calculations (with the pa- 
rameters obtained by fitting the low-energy bands only) 
are shown in Fig. 1121 The quantitative agreement is very 
satisfactory. 

The spin-orbit parameters entering the above fitting 
problem are in general momentum dependent, see gen- 
eralized k-dependent spin-orbit Hamiltonian in Table IVl 
of Appendix A. However, our numerical analysis shows 
that for the momentum vectors within the considered 
2.5% interval around the K point the values of the spin- 
orbit parameters modify less than 0.8%, the spin-orbit 



Hamiltonian at the K point, Eq. ([8]), is thus satisfactory. 

IV. CONCLUSIONS 

We have systematically investigated the spin-orbit cou- 
pling effects in bilayer graphene, both intrinsic and ex- 
trinsic, in the presence of a transverse electric field. We 
have presented first-principles results, and analyzed them 
from a tight-binding perspective. We have derived and 
used the most general spin-orbit Hamiltonian for the bi- 
layer graphene with an absent space inversion (gated bi- 
layer, bilayer on a substrate, or with adatoms) employing 
the underlying bilayer graphene symmetries. We have 
shown that the rough features of the spin splittings of 
the low-energy bands are well reproduced using a model 
with single-layer-like spin-orbit couplings and interlayer 
orbital hoppings. In particular, the intrinsic anticross- 
ings at zero external electric field are fully consistent with 
the anticrossing mechanism proposed for the single-layer 
graphene. In the presence of an electric field, the spin 
splittings of the (otherwise spin degenerate) bands are 
more subtle and complex than in a single layer graphene. 
Directly at the K point the value of the splitting of the 
low-energy bands is given by the intrinsic spin-orbit cou- 
pling. On the other hand, the high-energy bands are split 
in a proportion to the electric field, which is what is nor- 
mally expected. Away from the K points, our spin-orbit 
upgraded tight-binding model gives an excellent descrip- 
tion of the fine spin splittings due to electric field. With 
the help of the interlayer spin-orbit couplings we have 
fully covered the first-principles energetics at and near 
the K point. A quantitative and physically inspired fit- 
ting procedure was proposed to obtain the realistic val- 
ues of the interlayer spin-orbit coupling parameters; we 
have shown that these parameters have similar values as 
the intra-layer ones. The spin-orbit Hamiltonian we have 
proposed can be used in model studies of spin-dependent 
transport and spin relaxation in extended and confined 
bilayer graphene. 
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Appendix A: Effective spin-orbit Hamiltonian at 
K(K') from symmetry analysis 

The extrinsic spin-orbit coupling is induced by a trans- 
verse electric field which breaks the bilayer space inver- 
sion symmetry. Using symmetry group arguments we 
express the corresponding spin-orbit Hamiltonian at K 
in the basis of the on-site Bloch wave functions 4'Ai,s 
and ^'Bi,s [see Eqs. (j3|) and (jH])], in the format of Tab. 
im This Hamiltonian matrix has 10 real parameters 
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{•^11 7 -^n I ^^12 , , Ao = Ao + 2Abr, Aq = Aq — 2Abr, A4 = 
A4 + dXi, A4 = A4 — SX4, Ai, A3}. The first four parame- 
ters are related to the intrinsic intra-layer spin-orbit cou- 
phng Ai (hence subscript I); the rest (extrinsic) are la- 
beled corresponding to the hopping parameters 7^ with 
J = {0,1,3,4}. 



soc 


*Ai,t 


*Ai,4 


*Bi,t 


*Bi4 


*A2,t 


*A24 
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lAo 
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-lAo 
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lAo 












-Ai 


lAl 














AI2 







TABLE II: Spin-orbit matrix elements in bilayer graphene 
with external electric field at K. 



We now present the group theory analysis which leads 
to the spin-orbit coupling Hamiltonian in Tab. |lll In the 
absence of an electric field the point group symmetry of 
the bilayer unit cell is D^d, while the small group of K is 
D3. A transverse electric field along z breaks the space 
inversion symmetry. The point group reduces to C31, and 
the small group of K reduces to C3. The abelian group 
C3 = {i?, i?2ir/3j ^-27r/3} h^s three one-dimensional ir- 
reducible representations (we adopt the notation from 
Ref. Q), reproduced in Tab. IIIIl The complex con- 
jugate representations are related to the given ones by 
Ti = Ti and T2 = T3. 



representation 








r X r 




r2 




character xe 


1 


1 


1 












27r 


2n 




Ti 


r2 


Ts 


character xr 


1 


e' 3 




r2 












2w 












character Xr-^ 


1 


e ' 3 


e'— 




Ts 




r2 



TABLE III: Irreducible representations of C3 and their char- 
acters (left) and the corresponding multiplication table (right) 



The orbital components ^I'a; and ^E'Si of the on-site 
Bloch functions 4'Ai,s and ^'Si.s transform according to 
R2it/3i with the point-symmetry rotation axis given by 



A1-B2, as follows: 

R2^/3 (*A2 ) C-'^ ^-A^ , fl27r/3 (*B2 ) = *B2 ■ 

(AI) 

These transformations follow immediately from defini- 
tion, Eq. ([ij, when substituting for k- vector the K point 
momentum and employing the 7r-state property of 
orbitals. Hence 4'ai and belong to the trivial repre- 
sentation Ti, while '^Bi € r2 and ^'Aa to r3. Expressing 
the spin-orbit coupling in terms of and the spin raising 
and lowering matrices s± = ^{sx i iSy), 

Hso = -^^^{VVxp)-s - L,^s^+L-s++L+s- , (A2) 

we get operators and L± = L^^iLy which act only on 
the orbital part of the on-site Bloch wave functions 4'Ai.s 
and ^'Bi.s- The total potential V = Vbi + 21^1, entering 
Ec^. (|A2p , comprises the intrinsic bilayer potential energy 
and the energy due to the externally applied electric field. 
With respect to C3 the L operators transform as, 

-R27r/3-^2 = LzR2tt/3i R2tt/3L± = e^'~ -^±^^27r/3- (A3) 

Hence, € Fi, L_ G T2, and L_|_ G T^. The matrix 
element (vE'lLl'I)) belongs to the F.^ x Fl x F$ represen- 
tation of the group C3 . According to the (abelian) group 
theory: 

if Fvi, X Fl X F<{, 7^ Fi , than (*|L|$) = . 

Using this fact and the multiplication Tab. IIII[ one can 
identify all the vanishing Hamiltonian matrix elements in 
Tab. ini It is also clear that all diagonal matrix elements 
are real. We show below that the off-diagonal elements 
are either real or pure imaginary. 

The small group C3 of K is a subgroup of the 
point group C3V which in addition to the two rotations 
{i?-l-27r/3} also has three vertical reflections R±2Tr/3° 
Rv}- Here i?^ is the reflection with respect to the ver- 
tical plane defined by the electric field (z-axis) and the 
line connecting atoms Ai and Bi, see Fig. ([T]) for the co- 
ordinate system we use. Expressing the action of Ry in 
coordinates, 

Ry : r = (x, y, z) ^ {Ryv) = {~x, y, z) , 

we see that this transformation preserves the layer label 
and the pseudospin (A^ o A^ and B^ ^ Bi), but sends 
K to K'. The action of Ry on the on-site space Bloch 
wave functions and L operators is therefore as follows: 

-R«(«'a.,k) = «'A..K' = fl^K , RvLz = -L,Ry 

^«(*b.,k) = «'b.,K' = k > RvL± = L^Ry . 

(A4) 

Aloreover, 

(«'!$)= J d3r**(r)$(r) 

= J (i3r**(i?„r)$(i?,r) = (i?,„(^)|i?,($)). 
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that is, Rv is a unitary operator when acting on wave 
functions. The above facts suffice to show that, for ex- 
ample, (^fAi lizl'I'Ba) is real Indeed, 

Ai,K I i.*B..K> = ((i.*B„K)* I (*Ai,k)*> 

= (-i.(*B.,K)*|(*A,,K)*) 

= -(i.l'B^.K' |*Ai,K'> 

= -(i.i?„(*B.,K)|i?.(*Ai,K)> 

= -(-i?.„£.(*B..K) I i?.(«'Ai,K)) 

(i3*B2,K|*Ai,K) 
= (*Ai.K|i;^*B2,K)*- 

In a similar way one can show that the matrix elements, 
which comprise L± operators are imaginary. The only 
changes in the above computation appear in the second 
and the fourth lines. In the second line one uses {L±)* = 
—Lip and in the fourth line the commutation relation 
L::^Ry = RyL± instead LzRv = —RyLz- 

It is convenient to express the spin-orbit Hamiltonian 
given by Tab. |ll]in a more conventional form, 



HsO — Hi + HbR + -ffintcr + ^o: 



(A5) 



Here Hi and Hbr are the intrinsic and Bychkov-Rashba- 
like Hamiltonians, respectively. Hamiltonian iJintcr com- 
prises terms that are specific for interlayer coupling (bi- 
layer geometry), and H^i contains the remaining inter- 
layer contributions present only in a finite transverse elec- 
tric field. We introduce the Pauli matrices /x (layer), cr 
(sublattice pseudospin) , s (spin) , along with the unit 2x2 
matrices /zq, ctq, Sq, and symmetrized parameters: 



Ai = i(Ai2 - An) 
'^i = 5 (-^12 + An) 



K = 5(^12 + Ki) 



and 



2Abr, 
Ao 



5(^0 - Aq) 

^(Ao + Ag) 



A4 



i(A4 - A4) 

i(A4 + Al). 



We get the following effective Hamiltonian: 

i?i = i [(Ai - A[)cTa + (Ai + \[)az] HotSz 
+i [(Ai + A[)(To + (Ai - X'i)az] fi^rsz 



(A6) 



(A7) 



(A8) 



In the absence of an electric field, the spin-orbit Hamil- 
tonian of bilayer graphene was derived in Ref . [isj . This 
is a special case of our Hamiltonian in Tab. [Til We can 
restore space inversion symmetry 



I : r 



(x, y, z) i-> i(r) = (-a;, -y, -z), 



(A12) 



mapping K to K' and interchange the layer indexes and 
the pseudospin, 

Ai O B2 Bi o A2 . 

The orbital momentum L operators are invariant with 
respect to inversion, which is unitary, that is, | $) = 
(i(4') |i(<I')). The action of i on the on-site space Bloch 
wave function is: 



(A13) 



Applying the above facts to the nonzero matrix elements 
in Tab, [m we get the spin-orbit Hamiltonian of Ref. [l3| 
reducing the number of free (real) parameters to four, 



i(*Ai 


k) 


= *B2,K' 


= *B„K 


i(*B. 


k) 


= 4'Ai,K' 


= *li,K 


K*Bi 


k) 


= *A.,K' 


= *l2,K 


i(*A. 


k) 


= ^Bi.K' 


= *Bi,K 



{An = Aj;^,A 
given in Tab. 


12 — 

EYl 


X'12, Ao 




A4 = 


A4, Ai 


= A3 


= 0} 


SOC 




*Ai4 


*Bi.t 


*Bi4 


*A2,t 






*B2 4 


1'Ai,t 


A12 














1A4 








*Ai4 





— AI2 


lAo 

















*Bi,t 





— lAo 


-An 














— 1A4 


*Bi4 











An 














^A^.t 














An 











*A24 


— 1A4 














-An 


— iA() 






















lAo 


— A12 





*B24 








1A4 














AI2 



TABLE IV: Spin-orbit Hamiltonian matrix elements in bilayer 
graphene in the absence of an external electric field at K. 



Hbr - ^(AoMz + 2XBRfJ-o){Tcr^Sy - ays^) , (A9) example, that A3 = 0; by similar consid- 

erations one can check other spin-orbit elements: 



^^inter = -^{X^CTz + 5Xi(7o){T ^l^Sy + flyS^) , (AlO) 



(*Bi,K I i+*A3,K) = ((L+*A.,k)* I (*Bi,k)*) 

= (-i-(*A„K)*|(*B,,K)*) 
= -{L-^A2,K' |*Bi,K') 

Hcl = ^TSziUxCTx - ^J.yCry) ^^^^^ = - (i-i(^'Bi ,k) | i(^'A2,K)) 

+ ^[l^xira^Sy + aySa:) + tJ-y{T(7ySy - a^s^)] . = -(i-^-(*Bi.K) I K^'a^.k)) 

= -(L-'I'Bi.K I *A2,k) 

We have included t = 1 (K) and t = — 1 (K'), to cover ' ^ 

also K'. The translation of the Hamiltonian matrix ele- = — (^'Bi,k | L-*A2,k) 

mcnts from K to K' is based on time reversal symmetry. = — (5'bi,k | L+^a^,k) = . 
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Rewriting the matrix in Tab. IIVI in terms of the sym- 
metrized parameters as defined by Eqs. (|A6|) and (|A7p 
we get the spin-orbit couphng Hamiltonian in the ab- 
sence of electric field (that is in the presence of space 
inversion symmetry) in a more conventional notationi^^ 

Hso{E = 0) = AiT/ioo-^s^ + AiT/i^CTos^ 

A A (^^"^^ 

+ -ffJ-z{T<JxSy - (JyS^) - -fa-zirUxSy + HySx) ■ 

To explain quantitatively spin splittings of the low- 
energy conduction band, two spin-orbit parameters arc 
important: Ai and 4Abr = (Ao — Aq). The intcrlayer 
spin- flip parameters Ai, A3, A4 and SX4 arc important 
to quantitatively fit the valence band. Our spin-orbit 
Hamiltonian is valid in a general case of bilayer graphene 
with absent space inversion symmetry, such as gated bi- 
layer or bilayer on a substrate. It is likely that a stronger 
bonding substrate would make several of the ten param- 
eters much larger than they are in an applied field. 

The spin-orbit Hamiltonian Hso as presented in Ta- 
bic In] is, strictly speaking, valid directly at the K(K') 
point. It can be extended to be k-dependent. The strat- 
egy employs: the on-site Bloch states for a general k 
vector, Eqs. ([T]), the 7r-statc symmetry of p°J^ orbitals, 
transformational properties of the spin-orbit L operators. 



Eqs. (jA3|) . time-reversal symmetry, and finally, the near- 
est neighbor (nn) and the next-nearest neighbor (nnn) 
approximations, according to which: 

(vI/x,(k)|L|v|/Y^.(k)) « e"'(*-.-*'^-)x 

(A15) 

nn(n) 

In the expression above and Yj stand for an arbitrary 
couple of atoms Ai, Bi, A2, B2 in the bilayer elemen- 
tary cell. The underlying symmetries of bilayer graphene 
enable us to express: 

(pf (X,)|Lbf (Y,+R„„(„))) = c±'* {pfiX,)\L\pf{Y,)), 

where the phase factor is either 0, -I, or ^ depending on 
the atoms X^ and Yj. This observation simplifies sum- 
mation over Rjjn(jj) in Eq. (jAlSp . When proceeding as 
explained we arrive at the k-dependent spin-orbit Hamil- 
tonian shown in the Table FVl 

The structural spin-orbit functions u(k), v(k), w(k) 
and z(k) of the bilayer graphene are 



u(k) 
v(k) 
w(k) 
z(k) 



3^/3 



[sin k • Ri + sin k • R2 + sin k • R3] , 



3\/3 
1 



[sin k • Ri + e' 3 sin k • R2 -I- e~' 3 sin k • R3] 



[1 + e'' 



-+ikR2 



lpik-(tA2-tBi)-ik-R2 I'eik-Ra _,_ g- 



ikRs 



27r 



e' 3 



-ikRi 



(A16) 
(A17) 
(A18) 
(A19) 



where Ri = a(l,0), R2 = |-(— 1,— \/3) and R3 ~ ^{—1,^/3) arc Bravais hexagonal lattice vectors with the lattice 
constant a = 2.46A. 

The full tight-binding Hamiltonian with spin-orbit terms can be folded down to an effective Hamiltonian for the 
low-energy conduction and valence bands only. To this end we perform Lowdin transformatio n^ ''•'^^ by projecting 
Hamiltonian i?eff(k) given by Eq. ([9]), into the subspace of low-energy states with the on-site Bloch basis ^Si.t: 
^Bi4j ^A2,ti and ^A2.i- We keep only the intra-layer spin-orbit coupling parameters An, Aq and Aq as they are most 
relevant for the low-energy bands. The resulting effective 4x4 low-energy Hamiltonian, valid for k vectors close to 
the K point, for which 7o/(k) ^ 71 [the structural function /(k) is given by Eq. ([3])], is 



^eff^(k) 



/ V-Xn iAo^+/(k) 

-iAor;+r(k) V + Xn 
73/(k)* 

V-i(Ao/3+ + A[,r)/(k) 73/(k)* 



73/(k) 


-■^^ + Aii 
-iA[,rr/*(k) 



i(Ao/?+ + A^r)r(k) \ 
73/(k) 
iA^^-/(k) 
^V-Xn J 



(A20) 



with the dimensionless parameters 



/3^ = 



jV T A)7o ± 7174 

+ 72 - A2 
{V T A)74 ± 7o7i 



(A21) 



y2 



■7? 
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SOC(k) 


*Ai.t(k) 


*Ai4(k) 


*Bi,t(k) 


*Bi4(k) 


*A„t(k) 


*A,4(k) 


*B„t(k) 


*B24(k) 


*Ai,t(k) 


Ai2 u(k) 


iMAiv(k) 





iAq w (— kj 





1A4 w (k) 


Ai u(k) 


i^iv(k) 




-l^AiV (k) 


-Ai2u(k) 


\ /I \ 

lAo w(k) 





1A4 w(— k) 





-i^v (-k) 


A /I \ 

-Ai u(k) 


*Bi,t(k) 





— lAo w (k) 


-An u(k) 


-l^BiV(k) 





1A3Z (-k) 





— 1A4W (k) 


»T» /I \ 

*Bi,i(k) 


— lAo w(— k) 





l^BiV (k) 


An u(k) 


\ /IN 

1A3 z(k) 





— 1A4 w(— k) 





*A2,t(k) 





— 1A4W (— k) 





-1A3 z (k) 


An u(k) 


— l/iA2V(kj 





— iAqw (— kj 


*A,4(k) 


— 1A4 w(k) 





-iA3z(-k) 





lMA2V*(k) 


-Ki u(k) 


— lAg w(k) 





*B2,t(k) 


Aiu(k) 


i/iv(— k) 





1A4 w*(— k) 





iAq w* (k) 


-K2 u(k) 


i/iB2v(k) 


*B,4(k) 


— i/iv*(k) 


-Aiu(k) 


1A4 w(k) 





iAq w(— k) 





-l^B2V*(k) 


Ai2u(k) 



TABLE V: k-dependent spin-orbit matrix elements of the bilayer graphene Hamiltonian with external electric field within the 
nearest and next nearest neighbor approximation. Apart from the previously discussed spin-orbit parameters, there appear 
new five real couplings: fj,, hai, A*Bi, A*A2, and HB2, they are not present at the K(K') point. 
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